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Development of Methods of Predicting Large Crack 
Growth in Elastic-Plastic Work-hardening Materials 
in Fully Plastic Conditions 


Introduction and Objectives 

This report is in fulfilment of Research Grant Reference NAGAW - 3909 for NASA Ames 

Research Center, Dryden Flight Research Facility, Edwards, CA 93523. 

The objects of the first, exploratory, stage of the project were listed as:- 

f l.\ To make a detailed and critical review of the Boundary Element method as already 
published and with regard to elastic-plastic fracture mechanics, to assess its potential 
for handling present concepts in two-dimensional and three-dimensional cases. To 
this was subsequently added the Finite Volume method and certain aspects of the 
Finite Element method for comparative purposes. 

[ 2.} To assess the further steps needed to apply the methods so far developed to the 
general field, covering a practical range of geometries, work hardening materials and 
composites: to consider their application under higher temperature conditions. 

3. ; To re-assess the present stage of development of the energy dissipation rate, crack tip 

opening angle and J-integral models in relation to the possibilities of producing a 
unified technology with items 1 and 2. 

4. To report on the feasibility and promise of this combined approach and, if 
appropriate, make recommendations for the second stage aimed at developing a 
generalised crack growth technology for its application to real-life problems. 


Outline of the methods used 

To keep this exploratory study within reasonable bounds, only monotonic loading of metals 
at room temperature has been used in the case studies. Seven such studies have been made, 
four of the selected problems computed by all four methods used, followed by three problems 
specific to one or two of the selected methods. The four methods were ( 1) The Boundary 
Element [B.R]; (2) Finite Volume [F.V.]; (3) the Finite Element [F.E.] with a suite known 
as Modified WHAMSE [F.E(Wh)] and (4) the widely used ABAQUS program [F.E(Ab)]. A 
perhaps over-simplified statement on the differences in these approaches is that in the F.E 
formulation the classical differential relationships for equilibrium and compatibility are 
discretised for each nodal point; in the F.V. formulation these relationships are discretised in 
their integrated form for each sub-region in the domain (or control volume); in the B.E. 
formulation they are discretised in an integrated form for the whole domain. The B.E. and 
F.V. Methods are reviewed in Appendices A and B of the full report. 

All the programs are usable for linear elasticity and for incremental plasticity following the 
von-Mises criterion of yielding and the Prandtl-Reuss relationships, with work-hardening 
according to plastic potential theory, referred to hereafter as elastic-plastic. 
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The selected case studies 

1 . The cases-in-common to all methods used. 

The purpose of the four cases-in-common is to give some measure of the effort 
(preparation and computing time) and accuracy or suitability of results around which 
to make a discussion of the methods. A simply supported beam is used, with a central 
load, spread as a uniform pressure over a small distance on each side of the central 
point of contact. Two cases are elastic - (a) un-notched and (b) notched (Figs. 1 and 
2); and two elastic-plastic - (c) un-notched and (d) notched (Figs. 3 and 4). In all four 
cases plane strain is used. A coarse mesh "sighting shot" followed by medium and 
fine mesh computations were made to indicate convergence to about the 0.5% level. 

Sample results for the four cases are given in Tables 1 to 4, using fine mesh 
computations only together with a "datum” derived from other sources. There is no 
reason to suppose that the "datum" values are any more reliable than the calculations 
made by the team, but at least they give an idea of the correspondence with the best 
calculated and experimental data available to us. 

2. Further examples computed by particular methods are: 

a) With the F.V. and F.E(Wh) methods, the same beam as in Fig. 4 was 
computed for an impact loading at a constant striker speed of lOm/s. The results are 
given in the full report. 

b) Three-dimensional solutions of cases 1 to 4 were computed by 
F.E(WHAMSE), using a more refined mesh. The differences in values between the 
F.E(Wh) for the 3-D and 2-D cases are not large. The results are given in the full 
report 

c) Stable crack growth in an elastic-plastic beam. With the same beam as in 
Case 4, the stable growth of a crack was computed using F.E(Wh). Fig. 5 is a 
comparison of the computed and experimental load decay curves. 

The work is discussed in detail in the full report. 


Conclusions and recommendations for possible future work 

For elastic problems of monotonic loading in fracture mechanics, there is little doubt that 
linear elastic fracture mechanics is universally accepted as an appropriate model. Although 
fatigue is predominantly elastic in terms of the applied stress, it is now generally accepted 
that, for short cracks in relation to the micro-structure, it is inadequate, and inhomogeneities 
of micro-structure dominate the formation and early growth of a fatigue crack. 

For problems with plasticity there is less agreement but for many instances linear elastic 
fracture mechanics has been extended by the inclusion of a plastic zone correction factor. 

This is not the place to discuss the present stage of research into fracture, but merely to point 
out that there is a strong case for further computational methods by which the considerations 
of the physical models can be more easily carried on. 
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With the further development of these computational methods, it is expected that attention 
can be given to objectives 3 and 4 during the extension of the work recommended below. It 
will therefore be important to keep a "watching brier to incorporate such fracture studies at 
the right time. In particular, the J-integral procedure has been found wanting for both 
initiation and crack growth with extensive plasticity; initiation is being tackled by the Two 
Parameter methods, notably where crack tip fields are characterised by J and Q (hydrostatic 
stress); studies to describe crack growth by energy dissipation rate D, and associated crack 
tip opening angle (CTO A) exist, but are less well developed. We see a unification of J-Q and 
CTOA-D methods as the most promising approach and, as work proceeds, an additional 
project will be submitted to NASA when it is considered timely and promising to do so. 

In this regard, there is a considerable amount of useful work that can be done using the Finite 
Element method and there are many research groups pursuing this route. It is contended that, 
on the evidence of this preliminary study, there are clearly benefits to be realised from 
developing other methods and the Boundary Element and Finite Volume techniques are 
particularly promising. 

Our results show that the B.E method is an excellent choice for linear elastic stress analyses, 
including fracture problems. While some of its advantages over other methods are lost when 
non-linear material behaviour is introduced, it shows considerable promise for further 
development for elastic-plastic continuum (but not for composites). 

The Finite Volume method has potential for fracture problems, particularly non-linear such as 
in large plastic deformations or when they are dynamic. The method is also promising in 
regard to composite materials. 


In view of our preliminary findings, it is recommended that: 

1. Boundary Elements 

First, to explore the application of existing programs to two-dimensional elastic- 
plastic fracture mechanics problems. Attention would be given to the possible use of 
special crack-tip boundary elements to optimise the representation of the singular 
stress and strain fields there. There is also considerable scope for introducing 
iterative equation solution techniques of the congregate gradient type. The method 
should then be extended to three dimensional elastic-plastic problems. 

This would be a two-year project for two post doctorals working under Professor 
Fenner. 

2. Finite Volume method 

Future work on the F.V. method should concentrate on modelling non-linear 
constitutive relationships and the development of large strain formulation in fracture 
and should look at the extension to composite materials. 

A one to two year project is envisaged. 

3. Fracture Mechanics 

A continuing assessment of the progress of the computational work in order to start 
objectives 3 and 4 at the appropriate time. 


Table 1 

Computational Study No. 1 : Comparisons of fine mesh data only 
Un-cracked beam; elastic; plane strain 
Notation used (see Fig. 1) 
q/P compliance of C relative to D 

a(max) maximum bending stress at A (for a unit load) 


Results 

B.E. 

F.V. 

F-F./Whi 

F.E./ AM 

Datum 

q/P mm/MN 

31.09 

30.72 

31.17 

31.15 

31.03 [ 

a(max) MPa 

42.70 

42.52 

42.50 

42.72 

42.68* 

DoF 

272 

832 

1012 

1002 


(DoF ratios 

1 

3.0 

3.7 

3.7) 


CPU time, s 

0.3 

59 

100 

10 


(CPU ratios 

1 

180 

330 

30) 


CPU(s)/DoF 

0.0011 

0.071 

0.099 

0.010 



[1] Underwood etal, 1985. 

♦using conventional engineering elastic beam theory with the central load spread over 
±lmm as in the computations. 


Table 2 

Computational Study No. 2 : Comparisons of fine mesh data only 
Cracked beam; elastic; plane strain 
Notation used (see Fig. 2) 
q/P the compliance of C relative to D 
V/P the crack mouth opening compliance 
2v/P the compliances of the crack profile at two selected positions 
K[ the elastic crack tip intensity factor at unit load 


Results 

RJL 

F.V. 

F.E_.fWh) 

F.E./Abl 

Datum 

q/P mm/MN 

5.26 

5.28 

5.18 

5.24 

5.10 [1] 

V/P mm/MN 

3.24 

3.31 

3.17 

3.21 

3.26 [2] 

2/vP mm/MN 





327 [3] 

at a/32 

0.385 

0.433 

0.37 

0.381 

0.394 [2] 

at a/16 

0.563 

0.601 

0.53 

0.538 

0.568 [2] 

Kj MPa/ m 

0.936 

0.988 

0.964 

0.942 

0.953* 

DoF 

224 

1322 

6240 

6018 


(DoF ratios 

1 

5.5 

27.8 

26.9) 


CPU time, s 

1 

182 

900 

64 


(CPU ratios 

1 

182 

900 

64) 


CPU(s)/DoF 

0.0045 

0.15 

0.29 

0.011 



[1] Underwood et al [1985] 

[2] Gross etal [1968] 

[3] Kapp et al [1985] 

* using the lefm shape factor, Y = 10.650 from Srawley [1976], as listed by Towers [1981] in the 
form K = PY/B/W. 

The wide differences in degrees of freedom (DoF) and Central processor units (CPU) times between the 
methods are to be noted. They differ roughly by one to two orders of magnitude between the methods. The two 
F.E. methods (F.E.(Wh) and F.E.(Ab)) have very similar degrees of freedom but differ in CPU by near 15 fold 
for reasons of the method of solution used in this version of WHAMSE. The F.E.(Wh) data have been 
generated here by using an existing three-dimensional mesh run in plane strain rather than by generating an 
actual two-dimensional mesh. 


Table 3 

Computational Study No. 3 : Comparison of fine mesh data only 
Un-cracked beam; elastic-plastic; plane strain. 

Notation used below; (see Fig. 4.3) 

P(max) the final load at which a test run was stopped 
q displacement of C relative to D 
qpl plastic displacement of C relative to D 
q e l elastic displacement of C relative to D 

G(max) maximum bending stress at A. 


Results 

B.E. 

F.V. 

F.E(Wh) 

F.E.(Ab) 

Expt. [1] 

P(max) kN 

45.0 

45.0 

45.0 

45.0 

44.1* 

G(max) MPa 

1296 

1176 

1247 

1214 

1310 

q mm at P(max) 

4.95 

3.41 

5.45 

. 3.65 

10.0 

qpimm at P(max) 

3.55 

2.03 

4.05 

2.25 

8.1 

q e jmm at P(max) 

1.40 

1.38 

1.40 

1.40 

1.9 

qpjmm at 44kN 

3.47 

0.84 

2.71 

0.97 

8.0 

q D imm at 40kN 

0.95 

0.44 

0.55 

0.48 

2.7 

DoF 

344 

1100 

1908 

2314 


(DoF rados 

1 

3.2 

5.5 

6 . 7 ) 


CPU rime, m 

14 

130 

60 

27 


(CPU ratios 

1 

9.3 

4.2 

1.9) 


CPU Case 3/Case 1 

2800 

130 

35 

160 



[1] Dagbasi [1989]. * A loading roller slipped at this point whilst the load was still increasing. A 
maximum estimated by extrapolation is about 45 to 46kN; the test was of course in between plane 
stress and plane strain. 


Table 4 

Computational Study No. 4 : Comparisons of tine mesh data oniy 
Cracked beam; elastic-plastic; plane strain. 

Notation used (see Fig. 4.4):- 
q e l/P elastic compliance of C relative to D 
V e [/P elastic crack mouth opening compliance 
q. qpj displacement of C relative to D and plastic component thereof 
V, Vp| crack mouth opening displacement and plastic component thereof 


Results 

B.E. 

F.V. 

q e j/P mm/MN 

6.24 

6.25 

V e ]/P mm/MN 

3.93 

4.11 

Max q mm 

1.35 

1.37 

Maxqpimm 

0.02 

0.02 

Max V mm 

0.868 

0.905 

Max V p] mm 

0.019 

0.028 

F (max'' MN/m 

0.163 

0.133 

DoF 

264 

1944 

(DoF ratios 

1 

7 

CPU time, m 

5 

74 

(CPU ratios 

1 

15 

CPU(s)/DoF 

1.1 

2.3 

DoF Case 4/Case 2 

1.2 

1.5 

CPU Case 4/Case 2 

300 

25 


F.E(Wh) 

F.E.(Ab) 

Experiment [1] 
or Datum [2] [3] 

6.15 

6.20 

6.60 [1] 

6.63*; 6.03 [2] 

4.00 

4.06 

4.31 [1] 

4.47*; 4.07 [3] 

1.37 

1.37 

1.51 [1] 

0.04 

0.03 

0.03 [1] 

0.896 

0.904 

0.955 [1] 

0.032 

0.026 

0.024 [1] 

0.143 

0.136 

0.159 [1] 

6190 

4774 


23 

IB) 


50 

8 

0.5 

63 

13) 

0.8 


1.0 

0.7 


1.3 

60 



[1] Dagbasi, 1989. (Note the four computations here are made in plane strain whereas the beam 
tested is probably nearer plane stress; see above discussion and further results in Case 6. For the 
maximum load used in the computations, P = 216kN, a crack growth, Aa = 0.5mm is estimated from 
the experimental data whereas no growth was allowed in these studies). 

[2] Underwood et al, [1985], ^computed; plane stress and then plane strain. 

[3] Wu Shang-Xian, [1984], '‘computed; plane stress and then plane strain. 



Dimensions : 

S = 110 mm 
W = 13.52 mm 
B = 20.67 mm 
L = 120 &: 140 mm 


Properties : 

E = 2.0E+05 MPa 
Poisson’s ratio = 0.3 

Load. P = 1 kN 


r ig. 1 



Dimensions : 

S = 200 mm 
W = 50 mm 
B = 50 mm 
L = 250 mm 
a/W = 0.5 


Properties : 

E = 2.0E+05 MPa 
Poisson's ratio = 0. 

Load. P = 1 .0 kN 


Fig. ,2 
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Dimensions : 
W = 13.52 
S = 110 
= 120 
= 20.67 


L 

B 


w 


mm 

mm 

mm 

mm 
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Dimensions : 
a = 50.77 mm 

W = 94.2 mm 

S = 4W = 376.8 mm 

L = S+30 = 406.8 mm 
W B = 49.8 mm 
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